00100	TITLE FFT
00200	INTERNAL BITREV,FFT,FFTINT
00300
00400	BITREV:	0	;SUBROUTINE PUTS SPECTRA IN REVERSE BINARY ORDER
00500	BEGIN
00550	A←←0
00600	I←←1
00700	J←←2
00800	BITS←←3
00900	AC1←←3
01000	AC2←←4
01100		MOVEM 0,T0#
01200		MOVEM 1,T1#
01300		MOVEM 2,T2#
01400		MOVEM 3,T3#
01500		MOVE  0,0(16)
01600		HRRM C1
01700		HRRM C2
01800		HRRM C3
01900		HRRM C4
02000		HRRM C5
02100		HRRM  C6
02200		HRRM C7
02300		HRRM C8
02400		MOVE 0,@1(16)
02500		MOVEI I,1
02600		ROT I,@0
02700		MOVEM I,TWOMTH#	;2**M
02800		SOJ I,
02900		MOVEM I,J
03000	C:	CAMG I,J
03100		JRST B
03200		ROT I,1	;SWAP
03300		ROT J,1
03400	C1:	MOVE AC1,A(I)
03500	C2:	MOVE AC2,A(J)
03600	C3:	MOVEM AC2,A(I)
03700	C4:	MOVEM AC1,A(J)
03800		AOJ I,
03900		AOJ J,
04000	C5:	MOVE AC1,A(I)
04100	C6:	MOVE AC2,A(J)
04200	C7:	MOVEM AC2,A(I)
04300	C8:	MOVEM AC1,A(J)
04400		ASH I,-1
04500		ASH J,-1
04600	B:	MOVE BITS,TWOMTH	;DECREMENT
04700		ROT BITS,-1
04800		TDCN J,BITS
04900		JRST .-2
05000		SOJG I,C
05100		MOVE 0,T0
05200		MOVE 1,T1
05300		MOVE 2,T2
05400		MOVE 3,T3
05500		JRA 16,(16)
05600	BEND
05700
05800	;FAST FOURIER TRANSFORM
05900	FFT:	0
06000	BEGIN
06100	;SAVE THE ACCUMULATORS
06150		MOVE @(16)
06200		MOVEM IFS#
06300		MOVE [XWD 1,TEMPAC]
06400		BLT TEMPAC+17
06500		MOVE [XWD INITAC,1]
06600		BLT 17
06700
06800	;DEFINE ACCUMULATORS
06900	A←←0
07000	S←←0
07100
07200	AC1←←0
07300	AC2←←1
07400	T←AC3←←2
07500	TI←←AC4←←3
07600	UI←←4
07700	UR←←5
07800	I2←I3←←6
07900	;7
08000
08100	I←I1←←10
08200	N2←11
08300	J←FN←←12
08400	ILAST←FN2←←13
08500	LEXP←14
08600	LEXP1←15
08700	L←16
08800	NPL←17
08900
09000	;FOURIER ANALYSIS
09100		SKIPL IFS
09200		JRST SCL1
09300		MOVEI I,1
09400
09500	FALOOP:	FMPM FN,A+1(I)
09600		AOJ I,
09700	A0:	FMPM FN2,A+1(I)
09800		AOJ I,
09900		CAMG I,N2
10000		JRST FALOOP
10100
10200	;SPECIAL CASE L=1 LOOP
10300	SCL1:	MOVEI I2,2
10400	SCL2:	MOVE AC1,A(I2)
10500	A1:	MOVE AC2,A+2(I2)
10600	A2:	MOVE AC3,A+1(I2)
10700	A3:	MOVE AC4,A+3(I2)
10800	A4:	FADM AC2,A(I2)
10900	A5:	FSBM AC1,A+2(I2)
11000	A6:	FADM AC4,A+1(I2)
11100	A7:	FSBM AC3,A+3(I2)
11200		ADDI I2,4
11300		CAMG I2,N2
11400		JRST SCL2
11500
11600		MOVE M
11700		SOJE RET
11800
11900	;START OF L LOOP
12000	;START OF J=0 SPECIAL LOOP
12100	LLOOP:	MOVEI I,2
12200	SJLOOP:	MOVE I2,I
12300		ADD I2,LEXP1
12400		ADD I2,LEXP1
12500	A10:	MOVE AC1,A(I)
12600	A11:	MOVE AC2,A(I2)
12700	A12:	FSBM AC1,A(I2)
12800	A13:	FADM AC2,A(I)
12900	A14:	MOVE AC1,A+1(I)
13000	A15:	MOVE AC2,A+1(I2)
13100	A16:	FSBM AC1,A+1(I2)
13200	A17:	FADM AC2,A+1(I)
13300		ADD I2,LEXP1
13400		ADD I,LEXP1
13500	A20:	MOVE AC1,A+1(I1)
13600	A21:	MOVN AC2,A+1(I3)
13700	A22:	MOVE AC3,A(I3)
13800	A23:	MOVE AC4,A(I1)
13900	A24:	FADM AC3,A+1(I1)
14000	A25:	FADM AC2,A(I1)
14100		FSB AC4,AC2
14200	A26:	MOVEM AC4,A(I3)
14300		FSB AC1,AC3
14400	A27:	MOVEM AC1,A+1(I3)
14450		SUB I,LEXP1
14500		ADD I,LEXP
14600		CAMG I,N2
14700		JRST SJLOOP
14800
14900		MOVE L
15000		SUBI 2
15100		JUMPLE S120
15200
15300		MOVEM NPL,JJ
15400
15500	;START OF J LOOP
15600		MOVEI J,4
15700		MOVEI ILAST,4
15800		ADD ILAST,N2
15900		SUB ILAST,LEXP
16001	JLOOP:	MOVE AC2,JJ
16003		MOVN AC3,AC2
16005		ADD AC3,NT
16300	S00:	MOVE UR,S(AC3)
16400	S01:	MOVE UI,S(AC2)
16500
16600	;START OF I LOOP
16700		MOVE I,J
16800	ILOOP:	MOVE I2,I
16900		ADD I2,LEXP1
17000		ADD I2,LEXP1
17100	A30:	MOVE AC1,A(I2)
17200		MOVE T,AC1
17300	A31:	MOVE AC2,A+1(I2)
17400		MOVE TI,AC2
17500		FMP AC1,UI
17600		FMP T,UR
17700		FMP AC2,UI
17800		FMP TI,UR
17900		FSB T,AC2
18000		FAD TI,AC1
18100
18200	A32:	MOVE AC1,A(I)
18300	A33:	MOVE AC2,A+1(I)
18400		FSB AC1,T
18500		FSB AC2,TI
18600	A34:	FADM T,A(I)
18700	A35:	FADM TI,A+1(I)
18800	A36:	MOVEM AC1,A(I2)
18900	A37:	MOVEM AC2,A+1(I2)
19000		ADD I2,LEXP1
19100		ADD I,LEXP1
19200
19300	A47:	MOVE TI,A(I3)
19400		MOVN T,TI
19500	A50:	MOVE AC1,A+1(I3)
19600		MOVE AC2,AC1
19700		FMP TI,UR
19800		FMP T,UI
19900		FMP AC1,UR
20000		FMP AC2,UI
20100		FSB TI,AC2
20200		FSB T,AC1
20300
20400	A40:	MOVE AC1,A(I1)
20500	A41:	MOVE AC2,A+1(I1)
20600		FSB AC1,T
20700		FSB AC2,TI
20800	A42:	FADM T,A(I1)
20900	A43:	FADM TI,A+1(I1)
21000	A44:	MOVEM AC1,A(I3)
21100	A45:	MOVEM AC2,A+1(I3)
21150		SUB I,LEXP1
21200
21300		ADD I,LEXP
21400		CAMG I,ILAST
21500		JRST ILOOP
21600	;END OF I LOOP
21700		ADDM NPL,JJ
21800		ADDI ILAST,2
21900		ADDI J,2
22000		CAMG J,LEXP1
22100		JRST JLOOP
22200	;END OF J LOOP
22300
22400	S120:	ASH LEXP1,1
22500		ASH LEXP,1
22600		ASH NPL,-1
22700		AOJ L,
22800		CAMG L,M
22900		JRST LLOOP
23000	;END OF L LOOP
23100
23200		SKIPL IFS
23300		JRST RET
23400		MOVE N2,NN2
23500		MOVEI I,2
23600	A46:	MOVNS A+1(I)
23700		ADDI I,2
23800		CAMG I,N2
23900		JRST A46
24000
24100	RET:	MOVE [XWD TEMPAC,1]
24200		BLT 17
24300		MOVE IFS
24400		JRA 16,(16)
24500
24600	;FFTINIT (A,S,M,NT,N,N2,FN,NPL)
24700	↑FFTINT:	0
24800		MOVE (16)
24900		HLLI
24950		SUBI 2
25000		ADDM FALOOP
25100		ADDM SCL2
25200	FOR @$ DUMMY←0,50
25300	{ADDM A$DUMMY
25400	⎇
25500	MOVE 1(16)
25550		SUBI 1
25600	HRRM S00
25700	HRRM S01
25800	MOVE @2(16)
25900	MOVEM M
26000	MOVE @3(16)
26100	MOVEM NT
26200	MOVE @4(16)
26300	MOVEM NNN
26400		MOVE @5(16)
26500		MOVEM NN2
26600		MOVE @6(16)
26700		MOVEM TAC12
26800		MOVNM TAC13
26900		MOVE @7(16)
27000		MOVEM TAC17
27100		JRA 16,(16)
27200
27300	TEMPAC:	0
27400	BLOCK 20
27500	INITAC:	0 ↔ 0 ↔ 0 ↔ 0 ↔ 0 ↔ 0 ↔ 0 ↔ 0
27600	NN2:	0;	N2 11
27700	TAC12:	0;	FN 12
27800	TAC13:	0;	FN2 13
27900		10;	LEXP
28000		2;	LEXP1 15
28050		2;	L 16
28100	TAC17:	0;	NPL 17
28200	M:	0
28300	JJ:	0
28400	NT:	0
28500	NNN:	0
28600	BEND
28700	END